Load all required libraries.

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v ggplot2 3.3.5     v purrr   0.3.4
## v tibble  3.1.6     v dplyr   1.0.8
## v tidyr   1.2.0     v stringr 1.4.0
## v readr   2.1.2     v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(plotly)
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
library(broom)

Read in raw data from RDS.

raw_data <- readRDS("./year2.RDS")

Make a few small modifications to names and data for visualizations.

final_data <- raw_data %>% mutate(log_copy_per_L = log10(mean_copy_num_L)) %>%
  rename(Facility = wrf) %>%
  mutate(Facility = recode(Facility, 
                           "NO" = "WRF A",
                           "MI" = "WRF B",
                           "CC" = "WRF C"))

Seperate the data by gene target to ease layering in the final plot

#make three data layers
only_positives <<- subset(final_data, (!is.na(final_data$Facility)))
only_n1 <- subset(only_positives, target == "N1")
only_n2 <- subset(only_positives, target == "N2")
only_background <<-final_data %>% 
  select(c(date, cases_cum_clarke, new_cases_clarke, X7_day_ave_clarke)) %>%
  group_by(date) %>% summarise_if(is.numeric, mean)

#specify fun colors
background_color <- "#7570B3"
seven_day_ave_color <- "#E6AB02"
marker_colors <- c("N1" = '#1B9E77',"N2" ='#D95F02')
#remove facilty C for now
#only_n1 <- only_n1[!(only_n1$Facility == "WRF C"),]
#only_n2 <- only_n2[!(only_n2$Facility == "WRF C"),]

only_n1 <- only_n1[!(only_n1$Facility == "WRF A" & only_n1$date == "2020-11-02"), ]
only_n2 <- only_n2[!(only_n2$Facility == "WRF A" & only_n2$date == "2020-11-02"), ]

Build the main plot

      #first layer is the background epidemic curve
        p1 <- only_background %>%
              plotly::plot_ly() %>%
              plotly::add_trace(x = ~date, y = ~new_cases_clarke, 
                                type = "bar", 
                                hoverinfo = "text",
                                text = ~paste('</br> Date: ', date,
                                                     '</br> Daily Cases: ', new_cases_clarke),
                                alpha = 0.5,
                                name = "Daily Reported Cases",
                                color = background_color,
                                colors = background_color,
                                showlegend = FALSE) %>%
            layout(yaxis = list(title = "Clarke County Daily Cases", showline=TRUE)) %>%
            layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
        
        #renders the main plot layer two as seven day moving average
        p1 <- p1 %>% plotly::add_trace(x = ~date, y = ~X7_day_ave_clarke, 
                             type = "scatter",
                             mode = "lines",
                             hoverinfo = "text",
                            text = ~paste('</br> Date: ', date,
                                                     '</br> Seven-Day Moving Average: ', X7_day_ave_clarke),
                             name = "Seven Day Moving Average Athens",
                             line = list(color = seven_day_ave_color),
                             showlegend = FALSE)
      

        
        #renders the main plot layer three as positive target hits
        
        p2 <- plotly::plot_ly() %>%
          plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
                                       type = "scatter",
                                       mode = "markers",
                                       hoverinfo = "text",
                                       text = ~paste('</br> Date: ', date,
                                                     '</br> Facility: ', Facility,
                                                     '</br> Target: ', target,
                                                     '</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
                                       data = only_n1,
                                       symbol = ~Facility,
                                       marker = list(color = '#1B9E77', size = 8, opacity = 0.65),
                                       showlegend = FALSE) %>%
          plotly::add_trace(x = ~date, y = ~mean_copy_num_L,
                                       type = "scatter",
                                       mode = "markers",
                                       hoverinfo = "text",
                                       text = ~paste('</br> Date: ', date,
                                                     '</br> Facility: ', Facility,
                                                     '</br> Target: ', target,
                                                     '</br> Copies/L: ', round(mean_copy_num_L, digits = 2)),
                                       data = only_n2,
                                       symbol = ~Facility,
                                       marker = list(color = '#D95F02', size = 8, opacity = 0.65),
                                       showlegend = FALSE) %>%
            layout(yaxis = list(title = "SARS CoV-2 Copies/L", 
                                 showline = TRUE,
                                 type = "log",
                                 dtick = 1,
                                 automargin = TRUE)) %>%
            layout(legend = list(orientation = "h", x = 0.2, y = -0.3))
        
        #adds the limit of detection dashed line
        p2 <- p2 %>% plotly::add_segments(x = as.Date("2021-06-30"), 
                                          xend = ~max(date + 10), 
                                          y = 3571.429, yend = 3571.429,
                                          opacity = 0.35,
                                          line = list(color = "black", dash = "dash")) %>%
          layout(annotations = list(x = as.Date("2021-06-30"), y = 3.8, xref = "x", yref = "y", 
                                    text = "Limit of Detection", showarrow = FALSE))

        

        p1
        p2

Combine the two main plot pieces as a subplot

#seperate n1 and n2 frames by site
#n1
wrf_a_only_n1 <- subset(only_n1, Facility == "WRF A")
wrf_b_only_n1 <- subset(only_n1, Facility == "WRF B")
wrf_c_only_n1 <- subset(only_n1, Facility == "WRF C")

#n2
wrf_a_only_n2 <- subset(only_n2, Facility == "WRF A")
wrf_b_only_n2 <- subset(only_n2, Facility == "WRF B")
wrf_c_only_n2 <- subset(only_n2, Facility == "WRF C")


#rejoin the old data frames then seperate in to averages for each plant. 
wrfa_both <- full_join(wrf_a_only_n1, wrf_a_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))
## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke",
## "X7_day_ave_clarke", "Facility", "collection_num", "target",
## "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "se_L", "mean_total_copies",
## "sd_total_copies", "lo_95", "up_95", "log_copy_per_L")
wrfb_both <- full_join(wrf_b_only_n1, wrf_b_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))
## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke",
## "X7_day_ave_clarke", "Facility", "collection_num", "target",
## "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "se_L", "mean_total_copies",
## "sd_total_copies", "lo_95", "up_95", "log_copy_per_L")
wrfc_both <- full_join(wrf_c_only_n1, wrf_c_only_n2)%>%
  select(c(date, mean_total_copies)) %>%
  group_by(date) %>%
  summarize_if(is.numeric, mean) %>%
  ungroup() %>%
  mutate(log_total_copies_both = log10(mean_total_copies))
## Joining, by = c("date", "new_cases_clarke", "cases_cum_clarke",
## "X7_day_ave_clarke", "Facility", "collection_num", "target",
## "mean_copy_num_uL_rxn", "mean_copy_num_L", "sd_L", "se_L", "mean_total_copies",
## "sd_total_copies", "lo_95", "up_95", "log_copy_per_L")
#get max date
maxdate <- max(wrfa_both$date)
mindate <- min(wrfa_both$date)

Build loess smoothing figures figures

This makes the individual plots

#**************************************WRF A PLOT**********************************************
#add trendlines 
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_botha <- ggplot(wrfa_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_botha<<-..y..), method = "loess", color = '#1B9E77', 
              span = 0.25, n = 456)
## Warning: Ignoring unknown aesthetics: outfit
#look at the fits to align dates and total observations
#both
extract_botha
## `geom_smooth()` using formula 'y ~ x'

fit_botha
##   [1] 11.38408 11.43876 11.49244 11.54515 11.59692 11.64777 11.69774 11.74683
##   [9] 11.79508 11.84251 11.88915 11.93503 11.98016 12.02457 12.06829 12.11126
##  [17] 12.15340 12.19470 12.23513 12.27469 12.31335 12.35111 12.38794 12.42383
##  [25] 12.45877 12.49273 12.52571 12.55769 12.58864 12.61881 12.64838 12.67728
##  [33] 12.70546 12.73286 12.75940 12.78504 12.80970 12.83332 12.85584 12.87721
##  [41] 12.89734 12.91620 12.93370 12.95054 12.96732 12.98388 13.00005 13.01566
##  [49] 13.03054 13.04452 13.05742 13.06909 13.07935 13.08803 13.09496 13.09998
##  [57] 13.10291 13.10322 13.10068 13.09562 13.08832 13.07910 13.06826 13.05610
##  [65] 13.04293 13.02904 13.01476 13.00037 12.98619 12.97252 12.95967 12.94481
##  [73] 12.92539 12.90211 12.87564 12.84670 12.81596 12.78412 12.75188 12.71991
##  [81] 12.68892 12.65959 12.63262 12.60870 12.58851 12.56904 12.54708 12.52311
##  [89] 12.49760 12.47102 12.44383 12.41650 12.38950 12.36330 12.33836 12.31516
##  [97] 12.29416 12.27583 12.26063 12.24752 12.23513 12.22344 12.21245 12.20216
## [105] 12.19254 12.18360 12.17532 12.16770 12.16072 12.15438 12.14868 12.14359
## [113] 12.13911 12.13640 12.13636 12.13860 12.14273 12.14839 12.15518 12.16272
## [121] 12.17063 12.17853 12.18603 12.19274 12.19829 12.20230 12.20437 12.20512
## [129] 12.20547 12.20554 12.20543 12.20528 12.20518 12.20526 12.20563 12.20639
## [137] 12.20767 12.20958 12.21223 12.21574 12.22021 12.22532 12.23068 12.23630
## [145] 12.24222 12.24846 12.25504 12.26200 12.26935 12.27713 12.28536 12.29407
## [153] 12.30327 12.31301 12.32329 12.33416 12.34563 12.35773 12.37049 12.38393
## [161] 12.39999 12.42011 12.44361 12.46977 12.49789 12.52728 12.55724 12.58705
## [169] 12.61602 12.64345 12.66864 12.69089 12.70949 12.72374 12.73687 12.75235
## [177] 12.76983 12.78897 12.80943 12.83088 12.85296 12.87535 12.89769 12.91965
## [185] 12.94089 12.96106 12.97983 12.99685 13.01179 13.02430 13.03404 13.04068
## [193] 13.04386 13.04326 13.03852 13.02965 13.01712 13.00133 12.98265 12.96149
## [201] 12.93822 12.91324 12.88694 12.85970 12.83192 12.80397 12.77626 12.74916
## [209] 12.72308 12.69350 12.65646 12.61311 12.56461 12.51212 12.45680 12.39980
## [217] 12.34230 12.28545 12.23040 12.17833 12.13038 12.08771 12.05150 12.01629
## [225] 11.97684 11.93456 11.89088 11.84722 11.80500 11.76566 11.73060 11.70126
## [233] 11.67459 11.64687 11.61851 11.58992 11.56151 11.53370 11.50689 11.48151
## [241] 11.45795 11.43663 11.41797 11.40237 11.39025 11.38170 11.37632 11.37380
## [249] 11.37382 11.37608 11.38025 11.38602 11.39307 11.40110 11.40978 11.41881
## [257] 11.42786 11.43662 11.44478 11.45203 11.46108 11.47444 11.49149 11.51159
## [265] 11.53411 11.55841 11.58385 11.60982 11.63567 11.66078 11.68450 11.70620
## [273] 11.72526 11.74104 11.75582 11.77213 11.78976 11.80846 11.82802 11.84820
## [281] 11.86877 11.88950 11.91016 11.93053 11.95037 11.96945 11.98755 12.00444
## [289] 12.02051 12.03636 12.05199 12.06743 12.08268 12.09777 12.11271 12.12752
## [297] 12.14221 12.15680 12.17131 12.18575 12.20014 12.21450 12.22868 12.24257
## [305] 12.25618 12.26953 12.28264 12.29553 12.30824 12.32076 12.33314 12.34538
## [313] 12.35752 12.36956 12.38154 12.39346 12.40537 12.41726 12.42918 12.44113
## [321] 12.45313 12.46593 12.48008 12.49535 12.51152 12.52836 12.54563 12.56313
## [329] 12.58060 12.59784 12.61460 12.63067 12.64581 12.65980 12.67241 12.68341
## [337] 12.69294 12.70141 12.70908 12.71617 12.72291 12.72955 12.73630 12.74341
## [345] 12.75026 12.75614 12.76117 12.76544 12.76905 12.77210 12.77470 12.77695
## [353] 12.77894 12.78078 12.78258 12.78442 12.78642 12.78868 12.79110 12.79349
## [361] 12.79582 12.79806 12.80016 12.80209 12.80382 12.80530 12.80650 12.80739
## [369] 12.80792 12.80806 12.80778 12.80703 12.80473 12.80006 12.79344 12.78528
## [377] 12.77599 12.76597 12.75563 12.74540 12.73566 12.72684 12.71935 12.71358
## [385] 12.70996 12.70888 12.70933 12.70998 12.71085 12.71192 12.71320 12.71469
## [393] 12.71638 12.71829 12.72039 12.72270 12.72521 12.72793 12.73085 12.73396
## [401] 12.73829 12.74457 12.75244 12.76151 12.77142 12.78178 12.79221 12.80234
## [409] 12.81180 12.82019 12.82716 12.83232 12.83529 12.83570 12.83474 12.83379
## [417] 12.83280 12.83171 12.83046 12.82898 12.82721 12.82510 12.82259 12.81960
## [425] 12.81608 12.81198 12.80722 12.80175 12.79580 12.78962 12.78319 12.77648
## [433] 12.76947 12.76213 12.75444 12.74638 12.73792 12.72903 12.71970 12.70988
## [441] 12.69957 12.68874 12.67736 12.66540 12.65285 12.63968 12.62585 12.61121
## [449] 12.59565 12.57924 12.56208 12.54424 12.52581 12.50687 12.48750 12.46778
#assign fits to a vector
both_trenda <- fit_botha

#extract y min and max for each
limits_botha <- ggplot_build(extract_botha)$data
## `geom_smooth()` using formula 'y ~ x'
limits_botha <- as.data.frame(limits_botha)
both_ymina <- limits_botha$ymin
both_ymaxa <- limits_botha$ymax

#reassign dataframes (just to be safe)
work_botha <- wrfa_both

#fill in missing dates to smooth fits
work_botha <- work_botha %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_botha <- work_botha$date

#create a new smooth dataframe to layer
smooth_frame_botha <- data.frame(date_vec_botha, both_trenda, both_ymina, both_ymaxa)
#WRF A
#plot smooth frames
p_wrf_a <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_botha, y = ~both_trenda,
                    data = smooth_frame_botha,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_botha,
                                  '</br> Median Log Copies: ', round(both_trenda, digits = 2)),
                    line = list(color = '#1B9E77', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_botha, ymin = ~both_ymina, ymax = ~both_ymaxa,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_botha, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxa, digits = 2),
                                  '</br> Min Log Copies: ', round(both_ymina, digits = 2)),
                    name = "",
                    fillcolor = '#1B9E77',
                    line = list(color = '#1B9E77')) %>%
                layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF A") %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfa_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#1B9E77', size = 6, opacity = 0.65))

p_wrf_a
save(p_wrf_a, file = "./site_objects/wrf_a_year2.rda")
#**************************************WRF B PLOT**********************************************
#add trendlines 
#extract data from geom_smooth
#both extract
# *********************************span 0.6***********************************
#*****************Must always update the n = TOTAL NUMBER OF DAYS*************************
extract_bothb <- ggplot(wrfb_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_bothb<<-..y..), method = "loess", color = '#D95F02', 
              span = 0.25, n = 456)
## Warning: Ignoring unknown aesthetics: outfit
#look at the fits to align dates and total observations
#both
extract_bothb
## `geom_smooth()` using formula 'y ~ x'

fit_bothb
##   [1] 10.85586 10.93333 11.00953 11.08442 11.15798 11.23018 11.30100 11.37040
##   [9] 11.43835 11.50484 11.56982 11.63327 11.69517 11.75549 11.81419 11.87134
##  [17] 11.92702 11.98125 12.03406 12.08545 12.13545 12.18408 12.23135 12.27729
##  [25] 12.32192 12.36525 12.40730 12.44809 12.48765 12.52570 12.56203 12.59670
##  [33] 12.62979 12.66136 12.69149 12.72026 12.74773 12.77398 12.79907 12.82308
##  [41] 12.84609 12.86816 12.88936 12.90871 12.92535 12.93954 12.95152 12.96157
##  [49] 12.96992 12.97685 12.98260 12.98743 12.99161 12.99537 12.99899 13.00272
##  [57] 13.00681 13.01011 13.01141 13.01092 13.00883 13.00536 13.00069 12.99503
##  [65] 12.98858 12.98154 12.97411 12.96649 12.95889 12.95150 12.94452 12.93750
##  [73] 12.92986 12.92160 12.91274 12.90327 12.89321 12.88257 12.87135 12.85956
##  [81] 12.84721 12.83430 12.82085 12.80686 12.79234 12.77628 12.75790 12.73752
##  [89] 12.71546 12.69206 12.66764 12.64251 12.61702 12.59148 12.56622 12.54157
##  [97] 12.51784 12.49537 12.47449 12.45192 12.42475 12.39382 12.35999 12.32410
## [105] 12.28700 12.24954 12.21257 12.17693 12.14348 12.11306 12.08651 12.06469
## [113] 12.04845 12.03499 12.02108 12.00696 11.99286 11.97901 11.96562 11.95294
## [121] 11.94118 11.93057 11.92135 11.91373 11.90795 11.90423 11.90281 11.90369
## [129] 11.90665 11.91149 11.91805 11.92613 11.93556 11.94617 11.95777 11.97019
## [137] 11.98325 11.99676 12.01055 12.02445 12.03826 12.05401 12.07353 12.09645
## [145] 12.12237 12.15090 12.18165 12.21422 12.24822 12.28326 12.31894 12.35488
## [153] 12.39068 12.42595 12.46029 12.49332 12.52464 12.55386 12.58059 12.60443
## [161] 12.62813 12.65436 12.68262 12.71241 12.74325 12.77462 12.80604 12.83701
## [169] 12.86703 12.89560 12.92223 12.94643 12.96769 12.98552 13.00212 13.01985
## [177] 13.03850 13.05784 13.07765 13.09770 13.11778 13.13765 13.15710 13.17589
## [185] 13.19381 13.21063 13.22613 13.24009 13.25227 13.26246 13.27044 13.27597
## [193] 13.27884 13.27882 13.27569 13.26953 13.26072 13.24948 13.23603 13.22059
## [201] 13.20338 13.18463 13.16456 13.14338 13.12132 13.09860 13.07544 13.05207
## [209] 13.02869 13.00232 12.97034 12.93357 12.89283 12.84893 12.80268 12.75490
## [217] 12.70641 12.65801 12.61054 12.56479 12.52160 12.48176 12.44610 12.40970
## [225] 12.36809 12.32283 12.27544 12.22746 12.18043 12.13588 12.09534 12.06036
## [233] 12.02663 11.98935 11.94930 11.90727 11.86404 11.82040 11.77713 11.73501
## [241] 11.69484 11.65738 11.62344 11.59379 11.56922 11.54794 11.52769 11.50851
## [249] 11.49045 11.47355 11.45786 11.44342 11.43028 11.41849 11.40810 11.39914
## [257] 11.39167 11.38574 11.38138 11.37864 11.37874 11.38253 11.38957 11.39942
## [265] 11.41162 11.42572 11.44128 11.45785 11.47498 11.49222 11.50913 11.52525
## [273] 11.54013 11.55334 11.56701 11.58337 11.60204 11.62269 11.64494 11.66845
## [281] 11.69285 11.71779 11.74291 11.76785 11.79226 11.81578 11.83806 11.85873
## [289] 11.88011 11.90445 11.93126 11.96006 11.99038 12.02174 12.05366 12.08565
## [297] 12.11725 12.14797 12.17733 12.20485 12.23006 12.25247 12.27322 12.29375
## [305] 12.31404 12.33409 12.35389 12.37341 12.39266 12.41162 12.43028 12.44862
## [313] 12.46665 12.48434 12.50168 12.51867 12.53530 12.55154 12.56739 12.58285
## [321] 12.59789 12.61255 12.62688 12.64087 12.65456 12.66794 12.68103 12.69385
## [329] 12.70641 12.71871 12.73078 12.74262 12.75425 12.76568 12.77691 12.78798
## [337] 12.79815 12.80696 12.81477 12.82196 12.82889 12.83594 12.84348 12.85188
## [345] 12.86042 12.86817 12.87525 12.88175 12.88776 12.89339 12.89873 12.90387
## [353] 12.90892 12.91397 12.91912 12.92446 12.93010 12.93613 12.94250 12.94907
## [361] 12.95577 12.96257 12.96941 12.97625 12.98303 12.98971 12.99624 13.00256
## [369] 13.00863 13.01439 13.01981 13.02482 13.02935 13.03337 13.03696 13.04020
## [377] 13.04316 13.04592 13.04854 13.05110 13.05367 13.05633 13.05915 13.06220
## [385] 13.06556 13.06930 13.07357 13.07841 13.08367 13.08926 13.09503 13.10088
## [393] 13.10667 13.11230 13.11763 13.12255 13.12693 13.13066 13.13361 13.13566
## [401] 13.13777 13.14084 13.14461 13.14885 13.15330 13.15773 13.16188 13.16551
## [409] 13.16838 13.17024 13.17085 13.16997 13.16733 13.16271 13.15694 13.15095
## [417] 13.14467 13.13805 13.13101 13.12349 13.11540 13.10669 13.09728 13.08711
## [425] 13.07611 13.06421 13.05133 13.03742 13.02268 13.00736 12.99146 12.97496
## [433] 12.95785 12.94011 12.92173 12.90269 12.88299 12.86260 12.84152 12.81974
## [441] 12.79723 12.77398 12.74998 12.72523 12.69969 12.67337 12.64624 12.61821
## [449] 12.58920 12.55926 12.52845 12.49681 12.46440 12.43125 12.39743 12.36298
#assign fits to a vector
both_trendb <- fit_bothb

#extract y min and max for each
limits_bothb <- ggplot_build(extract_bothb)$data
## `geom_smooth()` using formula 'y ~ x'
limits_bothb <- as.data.frame(limits_bothb)
both_yminb <- limits_bothb$ymin
both_ymaxb <- limits_bothb$ymax

#reassign dataframes (just to be safe)
work_bothb <- wrfb_both

#fill in missing dates to smooth fits
work_bothb <- work_bothb %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothb <- work_bothb$date

#create a new smooth dataframe to layer
smooth_frame_bothb <- data.frame(date_vec_bothb, both_trendb, both_yminb, both_ymaxb)
#WRF B
#plot smooth frames
p_wrf_b <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_bothb, y = ~both_trendb,
                    data = smooth_frame_bothb,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothb,
                                  '</br> Median Log Copies: ', round(both_trendb, digits = 2)),
                    line = list(color = '#D95F02', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothb, ymin = ~both_yminb, ymax = ~both_ymaxb,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothb, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxb, digits = 2),
                                  '</br> Min Log Copies: ', round(both_yminb, digits = 2)),
                    name = "",
                    fillcolor = '#D95F02',
                    line = list(color = '#D95F02')) %>%
                layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF B") %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfb_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#D95F02', size = 6, opacity = 0.65))

p_wrf_b
save(p_wrf_b, file = "./site_objects/wrf_b_year2.rda")

#**************************************WRF C PLOT********************************************** #add trendlines #extract data from geom_smooth # *********************************span 0.6*********************************** #*****************Must always update the n = TOTAL NUMBER OF DAYS*************************

extract_bothc <- ggplot(wrfc_both, aes(x = date, y = log_total_copies_both)) + 
  stat_smooth(aes(outfit=fit_bothc<<-..y..), method = "loess", color = '#E7298A', 
              span = 0.25, n = 456)
## Warning: Ignoring unknown aesthetics: outfit
#look at the fits to align dates and total observations
#both
extract_bothc
## `geom_smooth()` using formula 'y ~ x'

fit_bothc
##   [1] 10.86752 10.92649 10.98448 11.04148 11.09749 11.15252 11.20654 11.25956
##   [9] 11.31157 11.36258 11.41257 11.46154 11.50948 11.55640 11.60228 11.64716
##  [17] 11.69108 11.73402 11.77598 11.81695 11.85694 11.89593 11.93392 11.97091
##  [25] 12.00689 12.04186 12.07581 12.10874 12.14064 12.17142 12.20103 12.22950
##  [33] 12.25685 12.28311 12.30832 12.33251 12.35570 12.37794 12.39924 12.41964
##  [41] 12.43917 12.45787 12.47575 12.49268 12.50850 12.52325 12.53696 12.54967
##  [49] 12.56142 12.57223 12.58214 12.59118 12.59939 12.60681 12.61346 12.61938
##  [57] 12.62461 12.62885 12.63181 12.63358 12.63422 12.63382 12.63245 12.63018
##  [65] 12.62708 12.62323 12.61870 12.61357 12.60792 12.60180 12.59531 12.58675
##  [73] 12.57473 12.55976 12.54235 12.52301 12.50224 12.48055 12.45845 12.43646
##  [81] 12.41508 12.39481 12.37618 12.35968 12.34583 12.33254 12.31759 12.30126
##  [89] 12.28384 12.26561 12.24687 12.22788 12.20895 12.19036 12.17238 12.15531
##  [97] 12.13944 12.12504 12.11241 12.10071 12.08896 12.07723 12.06556 12.05401
## [105] 12.04264 12.03152 12.02068 12.01020 12.00013 11.99052 11.98144 11.97294
## [113] 11.96507 11.95811 11.95220 11.94721 11.94298 11.93938 11.93629 11.93354
## [121] 11.93102 11.92858 11.92607 11.92338 11.92034 11.91683 11.91271 11.90758
## [129] 11.90131 11.89415 11.88636 11.87820 11.86993 11.86182 11.85411 11.84707
## [137] 11.84096 11.83603 11.83255 11.83078 11.83097 11.83205 11.83284 11.83346
## [145] 11.83403 11.83466 11.83546 11.83655 11.83804 11.84004 11.84267 11.84605
## [153] 11.85029 11.85550 11.86179 11.86929 11.87810 11.88834 11.90012 11.91356
## [161] 11.93047 11.95211 11.97774 12.00662 12.03804 12.07125 12.10553 12.14014
## [169] 12.17435 12.20743 12.23865 12.26727 12.29257 12.31382 12.33483 12.35963
## [177] 12.38774 12.41866 12.45193 12.48707 12.52360 12.56104 12.59891 12.63674
## [185] 12.67404 12.71035 12.74517 12.77804 12.80846 12.83598 12.86010 12.88036
## [193] 12.89626 12.90734 12.91311 12.91439 12.91248 12.90770 12.90032 12.89064
## [201] 12.87896 12.86557 12.85076 12.83482 12.81804 12.80073 12.78317 12.76565
## [209] 12.74847 12.72726 12.69827 12.66270 12.62173 12.57658 12.52845 12.47852
## [217] 12.42800 12.37810 12.33000 12.28491 12.24403 12.20856 12.17969 12.15211
## [225] 12.12055 12.08626 12.05053 12.01460 11.97976 11.94726 11.91837 11.89437
## [233] 11.87326 11.85231 11.83160 11.81123 11.79129 11.77188 11.75308 11.73500
## [241] 11.71773 11.70135 11.68597 11.67168 11.65857 11.64694 11.63697 11.62852
## [249] 11.62147 11.61570 11.61108 11.60750 11.60482 11.60292 11.60168 11.60098
## [257] 11.60069 11.60068 11.60084 11.60104 11.60234 11.60572 11.61093 11.61772
## [265] 11.62583 11.63503 11.64505 11.65566 11.66659 11.67760 11.68845 11.69887
## [273] 11.70862 11.71745 11.72621 11.73583 11.74621 11.75721 11.76872 11.78062
## [281] 11.79280 11.80513 11.81750 11.82978 11.84186 11.85363 11.86495 11.87571
## [289] 11.88673 11.89879 11.91174 11.92545 11.93975 11.95451 11.96957 11.98480
## [297] 12.00004 12.01515 12.02997 12.04438 12.05820 12.07131 12.08467 12.09924
## [305] 12.11486 12.13136 12.14856 12.16630 12.18441 12.20272 12.22106 12.23926
## [313] 12.25716 12.27457 12.29134 12.30730 12.32227 12.33609 12.34859 12.35960
## [321] 12.36894 12.37678 12.38344 12.38907 12.39377 12.39767 12.40090 12.40358
## [329] 12.40584 12.40780 12.40958 12.41130 12.41310 12.41509 12.41741 12.42016
## [337] 12.42276 12.42465 12.42602 12.42704 12.42791 12.42881 12.42992 12.43143
## [345] 12.43278 12.43336 12.43329 12.43268 12.43162 12.43024 12.42865 12.42695
## [353] 12.42526 12.42368 12.42232 12.42130 12.42072 12.42070 12.42022 12.41834
## [361] 12.41530 12.41135 12.40672 12.40165 12.39638 12.39115 12.38620 12.38177
## [369] 12.37809 12.37541 12.37397 12.37400 12.37510 12.37667 12.37869 12.38111
## [377] 12.38390 12.38701 12.39041 12.39407 12.39794 12.40199 12.40618 12.41047
## [385] 12.41483 12.41921 12.42454 12.43159 12.44007 12.44968 12.46013 12.47114
## [393] 12.48242 12.49368 12.50462 12.51497 12.52443 12.53270 12.53952 12.54457
## [401] 12.54900 12.55404 12.55952 12.56529 12.57118 12.57704 12.58272 12.58804
## [409] 12.59286 12.59701 12.60033 12.60268 12.60388 12.60378 12.60289 12.60178
## [417] 12.60043 12.59880 12.59686 12.59459 12.59195 12.58891 12.58544 12.58150
## [425] 12.57708 12.57213 12.56663 12.56054 12.55396 12.54701 12.53969 12.53198
## [433] 12.52387 12.51536 12.50643 12.49709 12.48732 12.47712 12.46647 12.45537
## [441] 12.44381 12.43179 12.41929 12.40630 12.39282 12.37885 12.36437 12.34933
## [449] 12.33373 12.31758 12.30091 12.28373 12.26608 12.24798 12.22945 12.21051
#assign fits to a vector
both_trendc <- fit_bothc

#extract y min and max for each
limits_bothc <- ggplot_build(extract_bothc)$data
## `geom_smooth()` using formula 'y ~ x'
limits_bothc <- as.data.frame(limits_bothc)
both_yminc <- limits_bothc$ymin
both_ymaxc <- limits_bothc$ymax

#reassign dataframes (just to be safe)
work_bothc <- wrfc_both

#fill in missing dates to smooth fits
work_bothc <- work_bothc %>% complete(date = seq(min(date), max(date), by = "1 day"))
date_vec_bothc <- work_bothc$date

#create a new smooth dataframe to layer
smooth_frame_bothc <- data.frame(date_vec_bothc, both_trendc, both_yminc, both_ymaxc)
#WRF C
#plot smooth frames
p_wrf_c <- plotly::plot_ly() %>%
  plotly::add_lines(x = ~date_vec_bothc, y = ~both_trendc,
                    data = smooth_frame_bothc,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothc,
                                  '</br> Median Log Copies: ', round(both_trendc, digits = 2)),
                    line = list(color = '#E7298A', size = 8, opacity = 0.65),
                    showlegend = FALSE) %>%
     layout(xaxis = list(range = c(mindate - 7, maxdate + 7))) %>% #buffer here
plotly::add_ribbons(x ~date_vec_bothc, ymin = ~both_yminc, ymax = ~both_ymaxc,
                    showlegend = FALSE,
                    opacity = 0.25,
                    hoverinfo = "text",
                    text = ~paste('</br> Date: ', date_vec_bothc, #leaving in case we want to change
                                  '</br> Max Log Copies: ', round(both_ymaxc, digits = 2),
                                  '</br> Min Log Copies: ', round(both_yminc, digits = 2)),
                    name = "",
                    fillcolor = '#E7298A',
                    line = list(color = '#E7298A')) %>%
                layout(yaxis = list(title = "Total Log10 SARS CoV-2 Copies", 
                                 showline = TRUE,
                                 automargin = TRUE)) %>%
                layout(xaxis = list(title = "Date")) %>%
                layout(title = "WRF C") %>%
  plotly::add_markers(x = ~date, y = ~log_total_copies_both,
                      data = wrfc_both,
                       hoverinfo = "text",
                       showlegend = FALSE,
                       text = ~paste('</br> Date: ', date, 
                                     '</br> Actual Log Copies: ', round(log_total_copies_both, digits = 2)),
                       marker = list(color = '#E7298A', size = 6, opacity = 0.65))

p_wrf_c
save(p_wrf_c, file = "./site_objects/wrf_c_year2.rda")

keeping in case

#save(wrfa_both, file = "./plotly_objs/wrfa_both.rda")
#save(wrfb_both, file = "./plotly_objs/wrfb_both.rda")
#save(wrfc_both, file = "./plotly_objs/wrfc_both.rda")
#save(date_vec_botha, file = "./plotly_objs/date_vec_botha.rda")
#save(date_vec_bothb, file = "./plotly_objs/date_vec_bothb.rda")
#save(date_vec_bothc, file = "./plotly_objs/date_vec_bothc.rda")
#save(both_ymina, file = "./plotly_objs/both_ymina.rda")
#save(both_ymaxa, file = "./plotly_objs/both_ymaxa.rda")

#save(both_yminb, file = "./plotly_objs/both_yminb.rda")
#save(both_ymaxb, file = "./plotly_objs/both_ymaxb.rda")

#save(both_yminc, file = "./plotly_objs/both_yminc.rda")
#save(both_ymaxc, file = "./plotly_objs/both_ymaxc.rda")